# !/usr/bin/env python
# -*- coding: utf-8 -*-
"""
@Time        : 2020/12/15 15:58
@Author      : Albert Darren
@Contact     : 2563491540@qq.com
@File        : single_variable_nonlinear_equation.py
@Version     : Version 1.0.0
@Description : TODO 自己实现单变量非线性方程组f(x)=0常见数值解法
@Created By  : PyCharm
"""
from sympy.abc import x
from sympy import init_printing, evalf, exp, solve, Eq

init_printing(use_unicode=True)


def bisection(f_x, interval_start: (float, int), interval_end: (float, int), epsilon=5e-3):
    """
    二分法求解单变量非线性方程组f(x)=0在区间[a,b]上的准确到小数点后指定精度的近似解
    :param f_x: 单变量非线性函数f(x)
    :param interval_start: 求根区间始点
    :param interval_end: 求根区间终点
    :param epsilon: 小数点精度
    :return: 近似解x,二分迭代次数k
    """
    if f_x.subs(x, interval_start) * f_x.subs(x, interval_end) < 0:
        # 初始化二分次数，区间长度
        bisection_count = 0
        interval_length = abs(interval_start - interval_end)
        while 1:
            middle_point = (interval_start + interval_end) / 2
            bisection_count += 1
            value = f_x.subs(x, middle_point)
            if interval_length <= epsilon or value == 0:  # 达到预定精度epsilon或者等于精确解时结束循环
                return middle_point, bisection_count
            if value * f_x.subs(x, interval_start) < 0:
                interval_end = middle_point
            else:
                interval_start = middle_point
            interval_length = abs(interval_start - interval_end)
    else:
        raise Exception("function has no solution in the specified interval")


def newtown_tangent(f_x, x0, true_root, epsilon=5e-4):
    """
    牛顿法求解单变量非线性方程组f(x)=0在区间[a,b]上的准确到小数点后指定精度的近似解
    :param f_x: 单变量非线性函数f(x)
    :param x0: 迭代初始近似值
    :param true_root: 精确解
    :param epsilon: 小数点精度
    :return: 近似解x,二分迭代次数k
    """
    f = f_x.subs(x, x0)
    diff_f_x = f_x.diff()
    diff_f = diff_f_x.subs(x, x0)
    iteration_count = 0
    while 1:
        x_k = x0 - f / diff_f
        iteration_count += 1
        if abs(x_k - true_root) < epsilon:
            return x_k, iteration_count
        x0 = x_k
        f = f_x.subs(x, x0)
        diff_f = diff_f_x.subs(x, x0)


def reduce_newtown():
    pass


def newtown_downhill():
    pass


def secant(f_x, x0, x1, true_root, epsilon=5e-4):
    """
    弦截法求解单变量非线性方程组f(x)=0在区间[a,b]上的准确到小数点后指定精度的近似解
    :param f_x: 单变量非线性函数f(x)
    :param x0: 第一个迭代初始近似值
    :param x1: 第二个迭代初始近似值
    :param true_root: 精确解
    :param epsilon: 小数点精度
    :return: 近似解x,二分迭代次数k
    """
    f_x0 = f_x.subs(x, x0)
    f_x1 = f_x.subs(x, x1)
    difference_quotient = (f_x1 - f_x0) / (x1 - x0)
    iteration_count = 0
    while 1:
        x_k = x1 - f_x1 / difference_quotient
        iteration_count += 1
        if abs(x_k - true_root) < epsilon:
            return x_k, iteration_count
        x0 = x1
        x1 = x_k
        f_x0 = f_x.subs(x, x0)
        f_x1 = f_x.subs(x, x1)
        difference_quotient = (f_x1 - f_x0) / (x1 - x0)


def parabola():
    pass


if __name__ == '__main__':
    # 第一组测试函数，来源详见李庆扬数值分析第5版P214，e.g.2

    # f = x ** 3 - x - 1
    # appropriate_root1, iteration_count1 = bisection(f, 1.0, 1.5)
    # print("二分法求方程近似解x={},迭代次数={}".format(appropriate_root1,iteration_count1))

    # 测试区间中点恰好是精确解时的情况
    # f = x - 1
    # appropriate_root2, iteration_count2 = bisection(f, 0.0, 2.0)
    # print("测试区间中点恰好是精确解时,近似解x={},迭代次数={}".format(appropriate_root2, iteration_count2))
    # 第二组测试函数，来源详见李庆扬数值分析第5版P238，e.g.7

    # f = x ** 3 - 3 * x - 1
    # initial_value = 2
    # true_solution = 1.87938524
    # root1, count1 = newtown_tangent(f, initial_value, true_solution)
    # print("牛顿切线法近似解x={},迭代次数={}".format(root1.evalf(), count1))
    # initial_value0, initial_value1 = 2, 1.9
    # root2, count2 = secant(f, initial_value0, initial_value1, true_solution)
    # print("弦截法近似解x={},迭代次数={}".format(root2.evalf(), count2))
    # 第三组测试函数，来源详见李庆扬数值分析第5版P224-225

    # c = 115
    # f = x ** 2 - c
    # true_solution = pow(c, 0.5)
    # initial_value = 11
    # solution, count = newtown_tangent(f, initial_value, true_solution, epsilon=1e-6)
    # print("牛顿切线法近似解x={},迭代次数={}".format(solution.evalf(), count))
    # 第四组测试函数，来源详见李庆扬数值分析第5版P229 e.g.10

    f = x * exp(x) - 1
    root = solve(Eq(f, 0), x)[0]
    print("调用sympy函数库解x={}".format(root))
    initial_value0, initial_value1 = 0.5, 0.6
    root, count = secant(f, initial_value0, initial_value1, root)
    print("牛顿切线法近似解x={},迭代次数={}".format(root.evalf(), count))
